A1 = readmatrix ('LCP_laser_intensity_spectra.csv')
A2 = readmatrix ('RCP_laser_intensity_spectra.csv')
B = readmatrix ('PDA10A2_responsivity.xlsx')

x1 = A1(9:83,1)
y1 = A1(9:83,4)

x2 = A2(9:83,1)
y2 = A2(9:83,4)

x3 = B(30:70,3)
y3 = B(30:70,4)

z = (y1-y2)./(y1+y2)*100

%photodetector responsitivity function
a1 = -4.876e-17
a2 = 3.398e-13
a3 = -8.011e-10
a4 = 8.979e-07
a5 = -0.0005266
a6 = 0.1581
a7 = -19.21

for i = 1:75
    n(i) = a1*x1(i)^6+a2*x1(i)^5+a3*x1(i)^4+a4*x1(i)^3+a5*x1(i)^2+a6*x1(i)+a7
end
n = n'
y1 = y1./n
y2 = y2./n

hold on
%plot(x1,y1,LineWidth=1,Color=[1 0 0])
%plot(x2,y2,LineWidth=1,Color=[0 0 1])
plot(x2,z,LineWidth=1,Color=[0 0 1])
xlim([500,750])
xticks(500:50:750)
ylim([-2.6,1.2])
yticks([-2:1:1])
xlabel('Wavelength [nm]')
ylabel(['Contast [%]'])
set(gca,'XMinorTick','on','YMinorTick','on')
set(gca,'box','on','Xcolor',[0 0 0],'Ycolor',[0 0 0])
hold off
box on
set(gcf,'Position',[0, 0, 200, 100])